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THEORY OF THE CORRECTION OF CELESTIAL OBSERVATIONS 
MADE FOR SPACE NAVIGATION OR TRAINING 
By Burnett L. Gadeberg and Kenneth C. White 
Ames Research Center 


SUMMARY 


This report discusses the theory of the correction of celestial 
observations, made for space navigation, instrument calibration, or training; 
as well as new theories concerning window refraction and lunar geometrical 
librations . It presents, in an easily understood form, essentially all 
current information on the corrections for celestial observations, providing 
a sound basis from which the engineer can derive high-precision observations. 
The corrections associated with these observations have been listed and 
classified. Where the theory has been completed, equations and constants 
are provided, or indicated, by which the reduction may be effected; in some 
cases, the equations have been recast in a form convenient for machine com- 
putation. Where the theory or knowledge was incomplete, developments are 
provided, or the requirements to develop the theory are made known. 


INTRODUCTION 


One of the requirements of a space navigation system is that celestial 
observations (measurements) of some kind be provided periodically, from which 
the state vector of the vehicle may be updated or determined directly. 

The types of measurements that may be classified as celestial 
observations are limited by standard usage in astronomy to angular measure- 
ments between celestial bodies. The body on which the observer is situated is 
not included, except when consideration is given to sextant horizon- altitude 
measurements. In general, however, we will be considering celestial observa- 
tions as angular measurements between various combinations of stars and 
bodies of the solar system, all situated some distance from the observer. 

Four types of observations are considered in this discussion. They are 
the measurements between: two stars; a star and the limb of a body displaying 

an extended disk; a star and some permanent marking on the surface, commonly 
known as a landmark; or two limbs of the same body or two different bodies. 

Any of these types of observations may be made with one setting of a sextant 
or two settings of a theodolite. All four types of observations are of value. 
To provide position information for space navigation, an observation must 
include at least one body within the solar system; however, star-star obser- 
vations are very useful for instrument calibration and for the training of 
observers. Measurements made between a star and the limb of an extended disk 
are useful to provide position information. Observations made between a star 
and a landmark also provide position information, but they may be more 



accurate when the destination body is close by. Sights between the two limbs 
of the same body are useful in providing distance information from the body 
when its apparent diameter is large. 

Such measurements have associated errors, which may be classified as 
either random or systematic. Accurate utilization of the data requires that 
these errors be compensated as precisely as possible. In the case of random 
errors, this compensation can be accomplished only by probabilistic means 
such as a weighted least squares (ref. 1) or a sequential data-processing 
method in the form of a Kalman filter (ref. 2). Studies describing such opti- 
mum compensation for random errors ignore the systematic errors on the premise 
that they are known and, therefore, can be corrected. However, the proper 
correction of systematic errors is not trivial and requires detailed attention 
to the theory . 

There has been much previous work, largely by astronomers, on the 
systematic errors in celestial observations (refs. 3 and 4). This work, how- 
ever, has been directed to satisfying particular needs, and does not consti- 
tute a complete theory for the problems peculiar to space navigation. One 
purpose of this report is to develop such a theory insofar as it is possible 
at this time, and to present mathematical techniques suitable for application 
to specific space-navigation problems. In some cases, astronomers’ previous 
work can be used directly, or with only minor modifications. Other cases 
require additional analysis, provided herein. Where there is insufficient 
knowledge for the development of precise correction formulas, the requirements 
of the correction theory are outlined in general terms. 

For most cases in space (position) navigation, the angle desired as an 
input to the system is the "geometric" angle as seen by the observer, where 
"geometric" angle means the instantaneous angle between the centers of the 
masses of the two bodies unaffected by distance, the speed of light, or the 
velocity of the observer. Since all the bodies of the universe are in a 
state of dynamical motion, and the angles between them are constantly changing 
for any observer, a measure of the time of observation is an essential factor. 

Even under ideal circumstances, an observer rarely can measure the 
desired angle. His own velocity causes aberration of the light, the apparent 
position of the body differs from the true position and ephemeris position due 
to the finite time of transmission of the light, the light rays may have been 
refracted by an atmosphere or a window, or the point on the body that really 
concerns him, the center of mass, is unmarked. There are also systematic 
errors associated directly with the use of optical instruments in making the 
measurements, and similar errors introduced by the human operator of the 
instrument. This report provides definitive descriptions of these errors and 
outlines some of the techniques whereby corrections may be applied. 

The results of this work have been used to develop a series of computer 
programs for reducing sights made on the surface of the Earth, in a high- 
flying airplane, or in a near-Earth satellite. These have been developed for 
the IBM 7094 digital computer, in FORTRAN IV language, and they correct 
the sights for aberration, refraction, librations, parallax, or changes of 
epochs. They have been used extensively in connection with studies made at 
Ames Research Center to evaluate the sextant -observer combination (refs. 5-7). 
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ANALYSIS OF CELESTIAL OBSERVATION CORRECTIONS 


General Procedure for Processing an Observation 

When observations are made for navigation purposes, the usual procedure 
for processing the observation data is to form a "residual," which is the 
difference between the actually observed quantity and the estimate of this 
quantity based on the best estimate of the observer's state (i.e., his posi- 
tion and velocity in space) at the time of the observation. The residual is 
then used, with appropriate weighting, to improve the estimate of the 
observer's state. When deterministic corrections are to be applied to the 
observation, the procedure consists of the following steps: 

(1) Obtain the estimated state of the observer (assumed given or 
computable) . 

(2) Obtain reference data for the sighting from star catalog and 
ephemeris , or both, etc. 

(3) Using (1) and (2), and the appropriate correction formulas, compute 
the corrected estimate of the line-of-sight vector. 

(4) From (3) , estimate the observation and form the residual by 
subtracting this from the actual observation (corrected for instrument error) . 


Classification of Corrections 

The corrections to be applied to celestial observations may be classified 
as follows: 

(1) Corrections due to finite velocity of light (aberration) , 

(2) Corrections dependent on position of the observer (parallax) , 

(3) Corrections dependent on time of observation (proper motion, 
ephemeris time. Earth precession and nutation, and lunar and planetary 
librations) , 

(4) Corrections due to medium of transmission (refraction) , 

(5) Corrections for instrument read-out errors (index, arc, and filter) , 

and 


(6) Miscellaneous corrections (irradiance, limb effects, semidiameter, 
and personal observer effects) . 

Of the six classes of observations, the first two, aberration and 
parallax, plus the fifth, instrument error, have been well defined by the 
astronomers. The only improvement here is to cast the equation for stellar 
aberration in vector form. The majority of the effects in (3) have been well 
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documented, but this report gives, for the first time, the correct vector 
matrix procedure for computing the geometrical lunar librations . Refraction 
through the atmosphere, part of (4), has been well defined. Refraction 
through an arbitrary window is developed here. The elements of (6) are still 
only poorly understood, and we will comment on what is known about them. 

Table I lists, for all of these types of corrections, the maximum 
correction magnitude to be expected for the four principal types of celestial 
observations (star-star, star-limb and occulations, star- landmark, and limb- 
limb) . It must be strictly understood that the values given are applicable 
only as qualified. Under different circumstances, corrections could be of 
vastly different magnitude than indicated. For example, local parallax is 
quoted as 3400", which is an approximate maximum value for an observer on the 
Earth sighting the Moon. But for other positions of the observer, the 
parallax may vary from 0° to 180°. 

Table I also provides information on the circumstances that affect the 
requirement for corrections. The applicability of the various corrections 
depends on: (1) the type of observation taken, (2) the type of mission, or 

more properly the mission profile variables at the time, and (3) the accuracy 
requirements of the sight. The first of these conditions may be readily 
determined from table I. For example, for a star-star sight made from space, 
the table shows immediately that planetary aberration, local parallax, ephem- 
eris time, Earth orientation, lunar librations, irradiance, limb effects, and 
semidiameter may all be omitted from consideration. Conditions (2) and (3) 
then determine the applicability of the remaining corrections, none of which 
may be safely omitted unless their size is known to be bounded within values 
that may be considered negligible for the specific case in question. 

A detailed analysis of each of the corrections and methods for applying 
them are given in the following sections. 


Corrections Due to the Finite Velocity of Light (Aberration) 

Aberration is the error in the apparent direction of the line of sight of 
an illuminated body due to the vector combination of the velocities of light 
and the observer relative to the inertial coordinate system. The direction 
of the observed light is rotated in the direction of the observer's relative 
velocity vector. 

In the field of astronomy, aberration effects are classified in three 
subdivisions: stellar (also known as annual), diurnal, and planetary. 

Stellar aberration is that component of aberration due to the velocity of the 
Earth in its orbit about the Sun (the observer assumed to be on the Earth) , 
and diurnal aberration the component due to the differential velocity imparted 
to the observer by the rotation of the Earth on its axis. Planetary aberra- 
tion is the term applied to the aberration phenomenon that exists in the 
sighting of relatively near bodies whose motion must be taken into 
consideration along with the motion of the observer. 
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TABLE I.- MAXIMUM ANGULAR CORRECTIONS DUE TO VARIOUS 
SOURCES FOR CONDITIONS SPECIFIED 



Star -star 

Star- limb and 
occultation 

Star- landmark 

Limb-limb 

Aberration 

Stellar t 

20 V 4 
for © 

20 V 4 
for © 

2074 
for ® 

N 

Planetary t t 

i 

0'.'67 
for ® B 

0767 
for ® B 

N 

Parallax 

Local + 

N 

3400" 
for © B 

3400" 
for © D 

31" 

for © D 

Stellar 

0778 

0'.' 7 8 

0778 

i 

Proper motion 

10'.'25/yr 

10'.'25/yr 

10725/yr 

i 

Ephemeris time f 

N 

18 V 1 

for ® D 

18V 1 

for © D 

N 

Earth 

orientation 

Precession 

N 

S0'.'27/yr 

S0'.'27/yr 

N 

Nutation 

N 

N <7 

N <7 

N 

Lunar 

librations 

Geometrical f 

i 

See 

irregularities 

11275 
for © 

See 

irregularities 

Physical =j= 

I 

See 

irregularities 

0V15 
for © 

See 

irregularities 

Re fraction 

Atmospheric O 

1860” 
for © 

I860" 
for ffi 

I860" 
for © 

4" 

for © 

Window 

* 

* 

* 

* 

Instrument 

Index 

A 

A 

A 

A 

Arc 

A 

A 

A 

A 

Filter 

A 

A 

A 

A 

Irradiance 

i 

1 16" | 

m 

1 32" 1 

Limb 

effects 

Atmospheric 

I 

cu 

N 

m 

Irregularities f 

i 

9” 

for s 

I 

18" 
for D 



Semidiameter f 

I 

930” 
for © D 

I 

1860" 
for © 3) 

Personal equation C 

20" 

20" 

20" 

20" 


® Earth 

B Moon 

* Magnitude depends upon window - can be computed from theory provided in this report 

I I Unknown, needs more research 

t Depends upon observer's velocity 

f Depends upon observer's distance 

A Depends upon instrument 

N Negligible, less than O'.'l 

I Inapplicable 

O Horizon observation from the Earth 

Cl As much as 20 sec of arc has been recorded for the personal equation for a highly 

trained observer. For an unskilled observer, the personal equation may be much 
higher. 

<7 The maximum coefficient in the expansions for nutation is 9V21. This has a negligible 

effect when correcting for atmospheric refraction, but may require consideration if 
observations are made of a near-Earth satellite. 









For our purposes, diurnal aberration will be computed simultaneously with 
the stellar aberration by simply using the total heliocentric velocity vector 
of the observer to compute the required correction. Note that for an observer 
on some other body (e.g., in a spacecraft) we simply use, in like manner, the 
total inertial heliocentric velocity vector of the observer's platform. 

Equations for computing the effect of stellar aberration are developed in 
appendix A. Either of two equations may be used, depending on the exactness 
desired. The rigorous equation is 


where 

s observed unit vector of the star 
S true heliocentric unit vector of the star 

V heliocentric unit velocity vector of the observer 

V velocity of the observer 

c velocity of light 

If this is not convenient to use, due to the absolute magnitude 
appearing in the denominator, a very close approximation is given by 

S = [i - $ tv • s)]s ♦ (|) v 

Planetary aberration needs to be treated in a different manner. Whereas 
the star catalog gives the direction of the unit vector to the star as it 
would be seen by a heliocentric observer, that is, the "true" line of sight, 
a planetary ephemeris gives the geometrical direction of the body as it 
actually exists from a point at the center of the Earth (moving in its orbit 
about the Sun) . A correction must be made because of the time required for 
light to travel from the observed body to the observer, as well as a correc- 
tion for the aberration. The correction for planetary aberration may be made 
by simply reducing the time of observation by the light time and taking the 
position of the body from the ephemeris at this reduced time. Thus, if an 
observation had been made of a planet at 09^30 m 20 s , the ephemeris position of 
the planet at this time would be used to determine the time required for the 
light to travel from the planet to the observer. If this light time were 
10 s , the time of observation would then be reduced to 09^30 m 10 s . The ephem- 
eris position of the planet at this reduced time would then give the correct 
position of the planet as seen by the observer at the time of observation. A 
more complete discussion of the phenomenon is given in appendix A. 
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Corrections Dependent on Position of Observer (Parallax) 


Star catalog and ephemeris data give the line of sight to celestial 
bodies from certain reference points (e.g., center of the Sun for star data) , 
but generally the actual observer is at some other point. This displacement 
of the observer from the reference point results in a change in direction of 
the incident light rays from the observed body - the phenomenon called 
parallax - and requires a correction to be applied to the data. 

The parallax corrections may be divided in two classes: one representing 

the parallax of relatively close objects, such as those within the solar 
system, and another representing stellar parallax. 


Local parallax (within the solar system) may, under certain circumstances, 
become very large. This type of parallax, however, may be handled very simply 
by making a linear vector transformation and thus presents no real obstacle to 

the determination of apparent position 
MOON of any known object. As an example, 
from sketch (a) , it may be seen that 
the vector position of the Moon may be 
corrected for parallax by writing the 
vector equation 



2 ) 


r = M - R 


( 1 ) 


Stellar parallax is in all cases 
very small. The maximum that has been 
recorded to date is only 0.78 second of 
arc, and that is for a relatively 
inconspicuous star of the eleventh 
Sketch (a) magnitude known as Proxima Centauri. 

The parallaxes of only a relatively few 
stars have been measured, due to the very remote positions of the stars in 
general, and most of these are of the order of 0 . 1 second of arc or less. 

Thus, it is a simple matter to avoid the use of stars with sensible parallaxes. 
If it is necessary to use a star with a high parallax, a linear vector 
transformation will account for the effect, as it does in the case of local 
parallax. 


Corrections Dependent on Time of Observation 

For economy in presentation, star catalog and ephemeris data are given 
only for certain reference times, and in using these data one must apply 
corrections that depend on the difference between the reference time and the 
time of the observation. These corrections are described below. 

Proper motion .- Proper motion is defined as the yearly motion of a star 
normal to the line of sight. In most cases, proper motion is very small. The 
largest that has been recorded belongs to a tenth magnitude star and is 
10.25 seconds of arc per year. Proper motion may be avoided by the judicious 
selection of stars for navigation purposes, the use of a star catalog of 
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recent epoch; or, since the proper motions are cataloged for both right ascen- 
sion and declination, simple corrections, linear with time, may be added to 
the reference stellar position. Possible nonlinear components of the motion 
need not be considered, for no proper motions have been measured to date with 
sufficient accuracy to indicate such effects. 

Ephemeris time.- Ephemeris time is the smooth flowing measure of time 
implied by~ the dynamical physical equations of motion, and is based princi- 
pally on orbital motions of the planets. All other measures of time commonly 
used in the field of astronomy, for example, universal time, are associated 
with the rotation of the Earth, which, for various and unknown reasons, dis- 
plays some irregularities. Ephemeris time is computed, for our purposes, from 

ET = UT + AT (2) 


where 

ET ephemeris time 

UT universal time (or Greenwich mean time) 

AT time interval 

The estimated value of AT, for a given year may be found in reference 8 
(p. vii) . The values of AT are essentially unpredictable but are only 
slowly varying from year to year. Final values for the past are determined 
from a large number of observations of the Sun, Moon, and planets. For the 
year 1900, the value of AT was -4 seconds of time and, for the year 1967, it 
is estimated to be 37 seconds of time. 

Although universal time is the standard used for civil reckoning and most 
astronomical computations, attention is here directed to ephemeris time 
because most ephemerides of the bodies of the solar system use it as argument. 
For example, if the transformation from universal to ephemeris time is neg- 
lected, the position of the Moon, as seen from the Earth, will be in error by 
approximately 18 seconds of arc. If the sight has been made in the lunar 
orbital plane, the error will show up directly as an error in the measured 
angle . 

Earth precession and nutation .- The orientation of the Earth (i.e., the 
equatorial coordinate system) changes with time and therefore is not neces- 
sarily the same at an observation time as at the reference time of the star 
catalog or ephemeris . The correct orientation may be important when correc- 
tions are being made for parallax and atmospheric refraction for an Earth- 
based observer. Also, data may be used that are expressed in coordinate 
systems based on the equatorial system for different epochs. In either case, 
corrections must be applied that amount to a transformation of coordinate 
systems . 

The gyroscopic forces that produce the Earth orientation changes result 
in complex motions containing secular terms and periodic terms. The secular 
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terms are easily represented by polynomial expansions in time, These are 
called precession, which has a maximum motion of approximately 50 seconds of 
arc per year. The aggregate of the periodic terms, which are of both long and 
short periods, are called nutation. 


A simple matrix equation permits the precession transformation from an 
initial epoch of 1950.0 to a later epoch: 

Cv Y y Z, 


V ! 


(3) 


Y, 


where 


V' 

X x 

"Xy 


any vector in the original coordinate system 
1.00000000 - 0 . 00029697T 1 2 - 0.00000013T 3 
-0.02234988T - 0.00000676T 2 - 0.00000221T 3 
-0.00971711T + 0.00000207T 2 + 0.00000096T 3 






(4) 


Y y = 1.00000000 - 0.00024976T 2 - 0.00000015T 3 


= Zy 


0 . 00010859T 2 - 0 . 00000003T 3 


Z z = 1.00000000 - 0.00004721T 2 + 0.00000002T 3 


J 


and T is the time interval between the two epochs expressed in Julian 
centuries of 36,525 days. Note that this matrix is orthogonal and that the 
reverse process may be obtained by using the transpose. 

The effects of nutation may be incorporated into the computations by 
premultiplying the above matrix for precession by the matrix 


1 

Al/i COS £ 

Alp sin e 


-Aip cos e 
1 

Ae 


-Aip sin e 
-Ae 
1 


(5) 


where 

e true obliquity of the ecliptic 

Aip nutation in longitude 

Ae nutation in obliquity 

and 
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£ 


e 0 + Ac 


E 0 = 23°27'08V26 - 46'.'845T - 0'.'0059T 2 + 0’.'00181T 3 > (6) 

Epoch 1900 Jan 0.313 / 

The two quantities Aip and &e must be tabulated values or computed from 
series expressions in which there are 

23 long-period terms \ 

1 for Aip 

46 short-period terms) 

16 long-period terms | 

> for Ae 

24 short-period terms ) 

Under the circumstances, one should critically review the accuracy 
requirements before it is decided to include corrections for nutation, since 
the coefficient of the largest term in nutation is only 9.21 seconds of arc. 
Since this represents a distance of only about 900 feet on the surface of the 
Earth, any effect on the atmospheric refraction corrections would be negli- 
gible. Similarly, for observations between natural members of the solar sys- 
tem or the stars, the effect of this 900-foot displacement would be entirely 
negligible. However, when observing an artificial satellite close to the 
observer, the effect may be quite large, and then the inclusion of the effects 
of nutation in the corrections would be appropriate. Most of the coordinate 
transformations will be made in conjunction with refraction and parallax 
corrections for more distant bodies, so one may conclude that, in all 
probability, nutation may be neglected. 

For more detailed discussions of precessions and nutation, one should 
consult references 8, 9, and 10 for the equations and constants presently in 
use, and reference 4 for a more detailed discussion of the theory and the 
development of the original constants. 

Lunar and planetary librations .- If landmarks are to be used in the 
navigation process, their locations must be known. Consider the lunar land- 
marks first. Since the Moon rotates about its own axis in the same period 
that it revolves in its orbit about the Earth, on the average the same face 
is presented to an observer on the Earth at all times. For this reason, the 
practice in recording the locations of lunar landmarks has been to refer them 
to a coordinate system associated with the mean Earth-Moon line. The presenta- 
tion of the Moon's face to the observer on the Earth is not exact, however, 
for three reasons. First, the Moon may be considered to rotate at a constant 
angular velocity about its axis, but its orbital angular velocity about the 
Earth is variable because of the eccentricity of the lunar orbit. This, com- 
bined with the inclination of the lunar equator to the orbital plane, causes 
an apparent oscillation in two dimensions about a mean position. These 
oscillations are known as geometrical librations since they are strictly a 
function of the Earth-Moon geometry existing at the time of observation. 
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Second, the figure of the Moon is that of a triaxial ellipsoid, which permits 
periodic rotational forces to be applied by the Earth and the Sun. These 
forces result in a complex pendulum- type oscillation known as physical libra- 
tions. Third, an observer on the Earth is displaced in position, during a 
12-hour period, by the diameter of the Earth. The resulting parallactic effect 
is known as diurnal libration. However, since diurnal libration is strictly a 
parallactic effect, it can be taken care of by the local parallax corrections. 

The method of handling the libration correction is to first compute, 
separately for both the geometrical and the physical librations, the three 
selenocentric angles describing the motion, and add the corresponding angles 
together. The unit vector to the lunar landmark then may be transformed 
through the three angles and the right ascension and the declination of the 
lunar position into a unit vector expressed in the equatorial inertial 
coordinate system. 

The equations used to compute the geometrical and physical librations of 
the Moon are relatively complex. Rather than introduce them at this time, it 
is expedient to show the complete development of the equations for the geo- 
metrical librations and to display them along with the final equations for the 
physical librations (as taken from ref. 9) in appendix B. 

In principle, observations of planetary landmarks could require a 
reduction similar to that used for lunar librations, except that planetary 
landmarks locations generally would not be expressed in terms of an Earth- 
planet mean line-of-sight-coordinate system as are lunar landmarks. Rather, 
they would usually be given in a planetocentric coordinate system, the orien- 
tation of which, in the observer's coordinate system, would have to be known 
(or computed) at any specific observation time to determine the line-of-sight 
direction. However, situations could occur in connection with planetary 
landmarks sightings such that the condition of commensurabi lity between the 
period of rotation and revolution exists (e.g., a synchronous satellite). 

Then a librational-type treatment would be reasonable for line-of-sight 
motions relative to a mean, and corrections similar to those for lunar 
librations would be required. 


Corrections Due to Medium of Transmission (Refraction) 

Atmospheric refraction .- Equations developed by astronomers for this 
purpose, such as those given in reference 11, are obtained by regression 
analysis techniques whereby empirical functional relationships are established 
between observed refraction and atmospheric variables. A major difficulty 
with such equations is that they are cumbersome to use. There is also the 
theoretical objection that, in reality, refraction depends on the total state 
of that portion of the atmosphere through which the light rays travel. How- 
ever, it is not practical to observe the total atmosphere, nor is it necessary, 
since the corrections obtained from the empirical relationships between 
refraction and surface conditions fit the observed data well enough for most 
purposes . 
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For these reasons, a simplified approximation to the empirical 
relationships is believed adequate for general use. One simplified equation 
(ref. 12), developed by Professor Comstock, approximates the refraction to 
within 1 second of arc for all altitudes of the observed body from the zenith 
down to 15° above the horizon except for extreme states of the atmosphere. 
This equation is 


r 


983b 
460 + T 


tgz 


( 7 ) 


where 

r refraction, seconds of arc 

b barometric pressure, in. Hg 

T atmospheric temperature at the observer, °F 
z zenith distance of the observed body, degrees of arc 

The refraction given by the value of r in the above equation lies in 
the vertical plane and apparently moves the observed body toward the zenith. 
This is due to the fact that the observer is situated at the center of the 
visible portion of the atmosphere, which has the shape of a planoconvex lens. 
Since any ray entering the eye of the observer must necessarily lie in a plane 
passing through the axis of rotation of this symmetry, the resulting refracted 
ray must also lie in the plane of symmetry. Thus, when making corrections for 
refraction, the ray must be transformed into the observer's topocentric 
coordinate system to determine the zenith distance of the observed body. 

It should be noted that the observer's local altazimuth system is fixed 
in relation to the equatorial system of the Earth and therefore undergoes con- 
tinuous precession and nutation along with the Earth. If an inertial system, 
or epoch, is used representing the equatorial coordinate system as of the 
beginning of the current year, the effects of precession and nutation on the 
refraction will be very small. The gyroscopic effects of the Earth can only 
affect the position of the observer's local coordinate system by approximately 
50 seconds of arc in any one year, and since the refraction is a slowly 
varying function, except when observations are made close to the horizon, the 
effect will, in general, be negligible. 

Window refraction .- The refraction of light rays occurring at the window 
of a space vehicle is a field that has not received scientific attention 
heretofore; consequently, although window refraction is amenable to mathemat- 
ical analysis, it has not been formalized to the same extent as atmospheric 
refraction. In general, two problems exist with window refraction. First, 
the design of the spacecraft often dictates that the windows be of irregular 
outline and that the spacecraft be pressurized. The resulting irregular 
three-dimensional curve of the window does not lend itself to the normal geo- 
metrical optics analysis, which assumes optical surfaces of revolution. 
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Second, an iterative procedure is involved in the determination of the vector 
from the sighting instrument to the window and thence to the particular 
celestial body under observation. A method of handling these calculations is 
outlined in appendix C; it is necessary to complete some preliminary research 
on the window before the computations may be carried out. This preliminary 
work consists basically of an analysis and testing of the deflection of the 
window under a differential load such that the matrix of the deflection curve 
polynomial may be found. A knowledge of the differential pressure and the 
orientation of the spacecraft at the time of observation then permits the 
calculation of the refraction of the line of sight of any given celestial 
body which may have been observed. 

As the equations and procedures for the reduction of the effects of 
window refraction are lengthy, their development and use is displayed in 
appendix C, which also compares the present theory with thin lens theory 
through the computation of a numerical example. 


Corrections for Instrument Read-Out Errors 

The sextant is used to make angular measurements between two celestial 
bodies, as needed for space navigation. The sextants of the past, and prob- 
ably of the future, are light, easy to use, and amazingly accurate, for their 
size and weight, compared to other astronomical instruments. The best results 
are obtained when the instrument is carefully adjusted, calibrated, and 
placed in the hands of a careful, well-trained observer. 

After the instrument is fully adjusted, three calibrations, or error 
corrections, remain to be made: index, arc, and filter. The index error is 

the reading of the instrument when both lines of sight are brought to bear on 
the same distant object. If the reading is not zero, an index error exists 
that must be applied to any sight obtained with the instrument. The arc error 
is the error made when the markings were scribed on the arc during its manu- 
facture. It is determined by sighting various bodies of known angular dis- 
placement and comparing these values with those read from the arc of the 
instrument. Filter errors are introduced when filters, used to reduce the 
intensity of bright images, are introduced into one of the lines of sight. 

Most filters have a slight wedge angle that deviates any ray of light being 
transmitted through them. The error is determined by measuring the position 
of an object before and after the filter is placed in the line of sight. One 
of the fullest discussions on the use, calibration, and adjustment of the 
sextant is given in reference 13. 

For accurate work the sextant must be very carefully calibrated. The 
index error must be determined sufficiently frequently that the magnitude, 
drift, and constancy are known within the limits of the precision required. 
With careful attention to observing technique, determination of the observer's 
personal equation, and calibration of the instrument, observations may be made 
to a precision greater than the least count of the instrument. 
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Miscellaneous Corrections 


Irradiance.- Irradiance is a psychophysiological effect wherein the 
image” of a bright body appears to be larger than its true diameter. It is 
apparently caused by an irregular refraction and scattering of light within 
the eye, which thereby stimulates the retina outside the immediate area 
covered by the illuminated image. This well-known effect was once thought to 
be caused by imperfect optics, as indicated in reference 13 (vol. II, p. 90), 
which was revised and last published in 1891. Although the phenomenon is 
now known to be principally psychophysiological in origin, very little is 
known about it and no satisfactory method has been devised to account for it 
at the present time. 

In reference 14, published recently, the effect was partially studied for 
very bright extended disks. Here irradiance effects from 38 seconds of arc to 
over 2° were measured. Reference 9 mentions corrections made for irradiance 
varying, in different cases, from 1.5 to 2.5 seconds of arc. No reference or 
explanations are given for these corrections, and the writers suspect that 
they are empirical, introduced for the purpose of rectifying consistent errors 
of observation (which is entirely legitimate) . One approach to reduce irra- 
diance is to cut down the apparent intensity of the extended body with neutral 
filters. However, additional research is needed here to indicate the effects. 
On the other hand, unpublished data of filtered sextant observations of the 
Moon (reduced and corrected by the authors) display an irradiance effect of 
15 seconds of arc. Thus it is apparent that the effects of irradiance, even 
for filtered observations, are large enough to seriously impair the required 
accuracy of celestial observations made for space navigation, and that 
insufficient data are available to correct for this effect. 

Limb effects .- Limb effects are similar to irradiance in that they affect 
the apparent diameter of the observed body. There are two sources of these 
effects. The first is the effect of an atmosphere on the observed body, which 
may cause an apparent increase in the visible diameter of the body, such as 
the Earth, when viewed from space. Very little is actually known about the 
atmospheric haze effect, except that it seems to be variable, and probably 
is affected by the terrestrial weather conditions. 

Second, the irregularities in the outline of the limb, such as the Moon, 
may cause either an apparent increase or decrease in the diameter of the body 
due to local topological features. These irregularities, due to the presence 
of mountains and valleys situated on the apparent limb, can only be conveni- 
ently accounted for, at the present time, when observations are made of the 
Moon by an observer in the vicinity of the Earth. They have been charted by 
Dr. C. B. Watts of the U.S. Naval Observatory (see ref. 15) for use with 
occultation measurements. These charts are also discussed in reference 9. 

For any other position of the observer, or for any body other than the Moon, 
this correction may not be considered since sufficient information is 
presently not available. 

Semi diameter . - Since many of the observations made for space navigation 
or testing will be made to the limb of an extended, illuminated body, a 
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correction for the semidiameter of the body must be applied to represent the 
observation to the center of the apparent disk. This may be accomplished by 
one of two methods. First, if the angular measure of the semidiameter is 
known for the observer's position, this may be directly added to or subtracted 
from the measured angle depending on whether the observation was made to the 
near limb or the far limb. Or, if parallax computations must be made, this 
correction may be incorporated directly into the linear vector transforma- 
tion, which was discussed under the section on parallax corrections. 

In addition to the correction for semidiameter, a small correction of 
-0.6 second of arc should be applied to the celestial latitude of the Moon to 
account for the misalinement of the center of mass and the center of figure of 
the Moon. Further reference may be made to reference 9 (p . 212) for this 
effect. For many situations, this last correction may be ignored because of 
its small size. 

Personal equation .- One form of systematic error that has had but little 
thought given to it is the error known to astronomers as the "personal equa- 
tion." The personal equation of an observer is that correction (error) which 
must be applied to his observations, after all other corrections have been 
made, to produce the true quantity that should have been observed. It may be 
observed as the mean residual error remaining after all other known correc- 
tions have been made. Personal equations may be quite sizable and consistent. 
In some cases, they have amounted to as much as 20 seconds of arc for a 
highly trained, very careful observer. Although the equation may be consis- 
tent, it generally varies over a period of time. It has been found to depend 
on, among other things, the type of instrument used, its size and position, 
the position of the observer's body, his degree of fatigue, and the type and 
velocity of the observed body. Any personal equation applied to a set of 
observations should have been determined in the recent past. It is obvious 
that a personal equation cannot be found when observing an unknown quantity. 
Thus, the personal equation should be determined during the training and 
retraining periods, and should be further updated by information obtained 
during the course of a mission. Further discussion on this subject may be 
found in references 12 and 13. 

At this point it should also be remarked that the philosophy of the 
observer in making high precision observations is also extremely important. 

High precision observations require not only a highly accurate, well- 
calibrated instrument, but an observer of knowledge, training, and a philos- 
ophy that will not allow him to do anything but his careful best. Only with 
such an observer will it be possible to have confidence in his personal 
weighting of individual observations. Discussions on the philosophy of weight- 
ing individual observations may generally be found in texts considering the 
method of least squares (refs. 13 and 16). When a number of sophisticated 
observations are made, weighting should always be used because a careful 
observer always has a feeling for the relative merits of his own observations. 
If he does not, he cannot be considered an accomplished observer and should 
be trained and practiced with this in mind until he can. 
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PROCEDURES FOR APPLYING CORRECTIONS 


Order of Making Corrections 

The order of the application of the corrections is important. A 
judicious order will eliminate second-order effects, reduce the computations 
involved, and eliminate the need for iterations. The method of ordering found 
to give the least confusion may be called the method of "reverse ordering." 
Consider the path of the light ray as it travels from the illuminated body 
to the eye of the observer. Along the path the various phenomena under 
discussion are encountered. Note the order of encounter with these phenomena 
and then apply their corrections in reverse order. If a case arises where a 
decision cannot be made between the time priorities for two of the corrections, 
it is considered good practice to apply the largest correction first; with 
the above rule of application it is apparent that the instrument corrections, 
index, arc, and filter should always be applied to the instrument reading 
first. On the other hand, since they all affect the light ray almost 
simultaneously, they may be applied in any order. 

Consider now the application of corrections to a sight made from the 
surface of the Earth between a star and the limb of the Moon. Tracing the ray 
from the star to the eye of the observer, the first phenomenon encountered by 
the ray is that of aberration as the light approaches the Earth, which is 
traveling through space with a velocity of the order of 18 miles per second. 
After this the ray is affected by refraction as it traverses the atmosphere 
to the observer. It then encounters the effects introduced by the sextant. 

The corrections are applied in the reverse order. The sextant corrections are 
first applied to the sextant reading. It is then appropriate to apply the 
corrections for refraction, which depend on the altitude angle as seen by the 
observer. To make these computations, one must use the known inertial coor- 
dinates of the star, the latitude and longitude of the observer, and the 
local sidereal time to transform the position of the star into the observer's 
local altazimuth coordinate system. 

After the refraction corrections have been completed, the aberration 
corrections may be applied as indicated. The above corrections must be applied 
to both lines of sight when more than one body is observed. However, the 
corrections applied to the two lines of sight are not necessarily the same, 
since, as indicated in table I, corrections for planetary aberration, libra- 
tions, irradiance, limb effects, and semidiameter are not applicable to a line 
of sight to a star. 

A step-by-step outline follows of the recommended procedure for each of 
the corrections, in the order in which they appear in table I. 

Aberration.- (stellar) 

1 . Compute the velocity vector of the observer in any convenient 
coordinate system. 
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2. Transform from that convenient coordinate system to the equatorial 
system of the desired epoch. 

3. Compute the vector of the star as it is seen displaced due to 
aberration using equation (Al) or equation (A4) . 

Aberration . - (planetary) 

1. Determine the distance of the local body from the fundamental center, 
using the recorded ephemeris . 

2. Using the parallax correction (linear vector transformation) 
determine the distance between the observer and the local body. 

3. Using the velocity of light, compute the time of passage of the light 
from the local body to the observer (light time) . 

4. Reduce the time of observation by the light time. 

5. Using the reduced time of observation reenter the ephemeris and 
determine the apparent position of the local body. 

Note : Before the ephemeris of the body may be entered, the time of 

observation and the reduced time of observation must be converted to ephemeris 
time (see section on ephemeris time) . 

Parallax . - (earth-based) 

1. Compute the radius vector of the observer from the center of the 
Earth using the equation for an oblate spheroid. 

2. For the vector of the observer in the geocentric equatorial system of 

date . 

3. Precess the vector from the equatorial system of date to the equa- 
torial system of the desired epoch (eq. (3)). 

4. Make the linear vector transformation as required for parallax 
(eq. (1)). 

Note: The computation of the radius vector of the observer from the 

center of the Earth must be made using the geocentric latitude of the observer. 

The formation of the vector of the observer in geocentric equatorial 
system of date requires the computation of the local mean sidereal time. The 
use of mean sidereal time is compatible with the use of precession and the 
omission of nutation. 

Parallax . - (non-earth-based) 

1. Compute the position vector from the fundamental center in any 
convenient coordinae system. 
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2. Transform from the convenient coordinate system to the equatorial 
system of the desired epoch. 

3. Make the linear vector transformation as required for parallax. 

Ephemeris time . - 

1. From the American Ephemeris and Nautical Almanac, select the AT 
correction for the year of the observation date. 

2. Calculate ephemeris time from the universal time by equation (2): 

ET = UT + AT. 

3. Enter the ephemeris as required. 

Note : The AT correction is given in the ephemeris only for the year of 

the publication and the past. If the correction is required for a future year, 
it will have to be estimated as best as possible by extrapolation of the 
values given. 

Past values of the AT correction are given to 0.01 second of time, and 
present values are given only to the nearest second due to the fact that the 
correction can only be determined by observation. As a consequence, it is 
really only known after the fact. 

Earth orientation .- (precession) 

1. Compute and form the transformation matrix from the polynomial 
expansions in time for the elements (eq. (4)). 

2. Transform the vector, using the matrix of equation (3) either from 
1950.0 to date, or using the inverse, from date to 1950.0. (Both transforma- 
tions will have to be used if the fundamental epoch of the computations is 
some date other than 1950.0.) 

Note : The polynomials for the elements and the matrix have been set up 
for the transformation from the epoch 1950.0 to the time used as the 
independent variable in the polynomials. 

The reverse of this procedure can be used to reverse the direction of 
the precession (transposed matrix) . 

The notation 1950.0 represents the beginning of the Besselian year, which 
is within a fraction of a day of January 1. It actually starts when the mean 
longitude of the Sun reaches exactly 280°. 

The fundamental epoch of the precession matrix (eq. (3)) may be changed 
to any other epoch by suitably adjusting the coefficients of equation (4) . 
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Librations . - (geometrical and physical) 


1. For the time of observation, compute the three libration angles using 
the matrix equations developed in appendix B for geometrical librations 

(eqs . (B5) and (B8)). 

2. Compute and add the increments to the libration angles using the 
equations for physical librations (eq. (B9)). 

3. Form the transformation matrix with the three libration angles 
(eq. (BIO)). 

4. Transform the selenographic vector of the landmark to the equatorial 
system of date. 

5. Precess the vector from the equatorial system of date to the 
equatorial system of epoch (eqs. (4) and (3)). 

Refraction . - (atmospheric) 

1. Transform the star coordinates from the equatorial system of the 
desired epoch to the equatorial system of date using the precession matrix 
(eqs. (3) and (4)). 

2. Transform from the equatorial system of date to the altazimuth system 
of the observer. 

3. Compute and apply refraction to the altitude of star (eq. (7)). 

4. Transform from the altazimuth system of the observer to the equatorial 
system of date. 

5. Transform from the equatorial system of date back to the equatorial 
system of the desired epoch using the transpose of the precession matrix. 

Note : The transformation between the equatorial system of date and the 

altazimuth system of the observer must be made with the geodetic latitude of 
the observer. 

This transformation also requires the computation of the local mean 
sidereal time from an expansion in mean solar time, which is compatible with 
the use of precession and the omission of nutation. 

The neglect of nutation produces only very small errors in the computed 
refraction . 

Refraction . - (window) 

1. Transform the observed star unit vector from the equatorial inertial 
coordinates to the window coordinates using equation (Cl) . 

2. Assume any convenient value of the internal unit vector V. (Values 
of a, 8, and y in eq . (CIO).) 
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3. Compute the intersection of this vector with the plane of the 
undeflected window using equation (Cll). 

4. Set V} = 0, equation (C12) . 

5. Using equations (03) to (06), iterate to find the intersection of 
the ray with the deflected surface of the window. 

6. Compute the normal to the surface using equation (07). 

7. Compute the refracted ray using equation (08). 

8. Repeat steps 3 to 7 for each refracting surface. 

9. Using equation (09), determine the angle between the free space 
refracted ray and the star vector determined in step 1. 

10. If the angle determined in step 9 is too large, compute a new 
internal unit vector using equation (C20) , including the gain factor computed 
from equation (C21) . 

11. Iterate the entire process from steps 3 through 10, until the angle 
determined in step 9 is small enough. 

12. The final refracted ray is then coincident with the star ray and 
the proper initial internal ray vector has been determined. 

Instrument . - (index, arc, filter) 

1. The arc correction should be determined first from the calibration 
curve of the sextant. Since the arc correction is a function of the location 
on the arc at which the reading is made, the arc correction should always be 
taken from the calibration curve for the actual reading of the sextant at the 
time of observation, that is, the recorded angle. 

2. The index correction and the arc correction may then be lumped 
together and applied to the sextant reading. 

3. The filter correction, if applicable, may then be applied. 

Note : It has been assumed here that the sextant has been properly 

adjusted for coplanarity of the principal plane of the sextant, the telescope, 
and the normal vectors of the horizon and index mirrors. 

Irradiance . - 

1. At the present time it is impossible to correct for irradiation 
effects when only limited observations are available. 

Note: The best that one can do presently is to use as heavy a filter as 

possible and hope for the best. The only possible result that can be expected 
is reduced accuracy. 
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Limb effects . - (atmospheric, irregularities) 


1. Atmospheric limb effects are very poorly known at the present time. 

2. Irregularities of the limb are in a condition impossible to use at 
the present time, except for some of the irregularities near the limb of the 
Moon, as seen from the Earth. 

Semidiameter . - 

1. From the relative position vector of the observed body, determined as 
indicated under parallax corrections, and from the known linear radius of the 
body, compute the semidiameter of the body. 

2. Add or subtract the above computed semidiameter from the measured 
angle according to whether the measurement has been made to the far limb or 
the near limb . 

3. Or, if direct measurements have been made of the semi diameter , step 1 
may be omitted. 


CONCLUDING REMARKS 


This report has presented new theories concerning window refraction and 
lunar geometrical librations. Further, we have indicated or collected in one 
location, and in an easily understood form, essentially all available informa- 
tion on the corrections for celestial observations, so that the knowledge 
required by the engineer for obtaining high precision observations is readily 
available. It has been pointed out that a careful analysis should be made of 
the mission profile to determine which of the many corrections may be safely 
omitted; some of the corrections require considerable computations and may 
amount to only a few seconds of arc in themselves, and may have an even smaller 
effect on the corrected observed angle. 


Ames Research Center 

National Aeronautics and Space Administration 

Moffett Field, Calif., 94035, Feb. 17, 1969 
125-17-05-01-00-21 
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APPENDIX A 


DEVELOPMENT OF VECTOR EQUATIONS FOR THE ABERRATION OF LIGHT 


Stellar Aberration 


The aberration of a star is the apparent displacement of the unit vector 
from the observer to the star due to the combination of the velocity of the 
observer with the velocity of light. This combination is the simple linear 
vector summation of the two velocities, reduced to the unit vector form. 

Since the velocity of light in free space is considered a constant, if 
the velocity vector of the observer were also constant, all the apparent posi- 
tions of the stars would be displaced by different but constant angles from 
their true positions. Under these conditions a slightly erroneous but con- 
stant picture of the heavens would be obtained. Now the velocity of the Sun, 
or the solar system, is essentially constant, at least over very long periods 
of time, within the accuracy of our astronomical measurements. Consequently, 
all star catalogs quote star positions as though they were seen by an helio- 
centric observer moving with the velocity of the Sun. The important point to 
be remembered here is that the correction for aberration must be based on the 
heliocentric velocity of the observer. 


Two vector equations will be developed for the aberration of light, one 
will be rigorous, and one will be slightly more convenient to use for machine 
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Figure 1.- The aberration of light due to the 
heliocentric velocity of the observer. 


computations. The geometry of the 
situation is displayed in figure 1. An 
observer, located at 0 and moving 
through space with the velocity vector 
V, observes a star whose true direction 
(heliocentric direction) is given by 
the unit vector S. The light from the 
star approaches the observer with the 
velocity vector c" relative to fixed 
space. The velocity of the star light 
relative to the moving observer is 
given by -s~, and the star appears to 
the observer to be located in direction 
s. Now, if we use the convention that 
any vector is composed of a magnitude 
and a unit vector, then it is apparent 
that 

V = vv 

c = cc = -cS 
Then we may write 

s = V - c 
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and 


VV + cS 


or 


V - c 

s = _ 

| V - c | | VV + cS | 

s = CV/c)V + S 
(V/c)V + S 


(Al) 


This is the rigorous equation for stellar aberration. 1 It has the disadvantage 
that the magnitude of the vector must be evaluated for the denominator. A 
slightly more convenient equation may be derived from figure 1 with the help 
of two additional vectors. These auxiliary vectors are a" and x, both of 
which are normal to the direction to the star as given by S. Then w may 
write 

s = S + a (A2) 

now 

a << 1 

since 


x 

a = — 
c 


(A3) 


and 


V • S << c 


Then we may also write 

x = S x (V x S) 


which may be reduced to 


x = -V(V • S)S + VV 


^hese equations have been developed in nonrelativistic mechanics and 
are thus accurate for velocities of the observer where V/c << 1. It has been 
kindly pointed out by Mr. William M. Adams, Jr., of Langley Research Center, 
that the correct equation for a relativistic observer would be 



This equation reduces to equation (A4) for nonrelativistic observers where 
(V/c) 2 and higher terms may be neglected. 
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Then, substituting into equation (A3) , 

a = - - (V • S)S + - V 
c c 


and then into equation (A2) 

s 


1 


(V • S) 



(A4) 


Although equation (A4) is an approximation, it is a very good approxima- 
tion, since the values of V are always very small compared to the velocity 
of light c. The value of tg _1 (V/c) for the Earth, which is known as the 

constant of aberration, is only 20.4 seconds of arc. 


Planetary Aberration 

When an observer views a planet, or any other body within the solar 
system, he sees the body displaced by stellar aberration just as he would 
when viewing a star. However, a fundamental difference exists in the method 
of correcting for this error (see ref. 17). This difference arises as a 
direct result of the fact that the position of the planet, as supplied by the 
ephemeris, is fundamentally different information than that supplied by a star 
catalog for the position of the star. A star catalog gives the direction from 
which the light arrives when it enters the eye of the heliocentric observer. 

The fact that the observer is considered to be heliocentric is not really 
germane to this immediate discussion, but the fact that the catalog supplies 
the direction from which the light is approaching is important. The basic 
ephemeris of a planet, on the other hand, gives the true geometric position of 
the planet in space at the time specified, and makes no allowances for the 

travel time of the light or 
any other phenomenon. Although 
this may give the impression 
that some complication would 
be involved in the computation 
of the correction, it will be 
shown that the apparent posi- 
tion of the planet may be 
taken directly from the ephem- 
eris by the simple expedient 
of decrementing the time of 
observation by the travel time 
of the light from the planet 
to the observer. 

OBSERVER^ Figure 2 shows the 

PATH combined effects of both 

stellar and planetary aberra- 
tion. At the time the light 
leaves the body at M 0 (say, 

Figure 2.- Planetary and stellar aberration combined. the Moon) the observer is 

located 0 o . When the light 
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reaches the observer at 0^, along the path r, the 
Moon has moved on to . If at the time of observa- 
tion we entered the ephemeris of the Moon it would 
indicate that the Moon was located at Mj at the 
distance r m from the observer. Now, the stellar 
aberration discussed above gives the observer the 
impression that the Moon is located in the direction 
of the unit vector s. We assume that during the 
time of light passage from M 0 to 0^ that the 
observer has moved with uniform velocity along d 
from 0 o to 0 X and arrives 0^ with the velocity 
vector V, and that the approaching light has the 
velocity vector c\ Then, as discussed above in 
conjunction with equation (Al) , the body appears to 
be located in the direction s, which is parallel to 
s". Now it is apparent that "c and T are colinear, 
as are V and d. Then we may write, using At for 
the light time 


then, on elimination of At, 

V = d 
c r 


Thus, the two triangles are similar and the side D 
is parallel to the side s'. Now s or s is the 
apparent direction of the Moon, as seen by the 
observer, and D is the direction of the Moon from 
the observer at the time the light left the Moon. 
Consequently, the apparent direction of the Moon, at 
the time of observation, is the same as the true 
direction of the Moon at the time the light left the 
Moon. At the time of observation, the light has 
traveled the distance r in the time At. However, 
neither of these quantities is known. We do know the 
distance r m from the ephemeris, and we may assume 
that 



since the velocity of the Moon is very small compared 
to the velocity of light and consequently r and r m 
are approximately equal. Thus, the apparent direc- 
tion of the Moon at the time of observation may be 
found by entering the ephemeris at the time of obser- 
vation, determining the distance to the Moon and the 
light time, and then reducing the time of observation 
by the light time and then reentering the ephemeris 
to determine the apparent direction of the Moon at 
the time of observation. 
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APPENDIX B 


DEVELOPMENT OF EQUATIONS FOR LUNAR GEOMETRICAL LI BRAT IONS 


The rotation of the Moon about its axis may be separated into two 
components. The first of these is a steady rotation in a counterclockwise 
direction, as seen from the north, which is in the same direction as the 
revolution of the Moon about the Earth in its orbit. The period of this rota- 
tion is exactly equal to the sidereal period of the lunar orbit. This rota- 
tion, combined with the orbital motion, the inclination of the Moon's equator 
to the ecliptic, the position of the lunar node, and the inclination of the 
ecliptic to the Earth's equatorial plane gives rise to the geometrical libra- 
tions . The second part of the lunar rotation is a very low amplitude oscil- 
latory motion, which gives rise to the physical librations . It is the 
geometrical librations that are of concern here. 

Several pertinent facts concerning the interrelationship between the 
rotation and revolution of the Moon are expressed in the three laws of Cassini. 
These three laws, which were originally empirical and later proved 
mathematically by Lagrange and Laplace, state that: 

1. The Moon rotates eastward, on an axis fixed in its body, with a 
constant angular velocity, and the period of rotation is exactly the same as 
that of the lunar sidereal revolution about the Earth. 

2. The inclination of the lunar equator to the ecliptic is a constant. 

3. The poles of the lunar equator, the ecliptic, and the lunar orbit all 
lie in a great circle, in that order. 

The rotation of the Moon about its polar axis is constant and the angle 
through which it has turned after a given period of time is a linear function 
of the time. Consequently, there is a direct relationship with the mean 
longitude of the lunar orbital motion. If the mean motion of the Moon and the 
true motion coincided, then exactly the same face of the Moon would be pre- 
sented toward the Earth at all times. However, the osculating orbit of the 
Moon is eccentric, and it is also perturbed by the Sun and the planets. As a 
consequence, the difference between the true position and the mean position 
causes the Moon to appear to oscillate about the radius vector of the Moon 
from the Earth. 

The investigation will now proceed along two parallel paths. In one we 
will develop the transformation matrix for any s elenographic vector to the 
equatorial system via the laws of Cassini and the general geometry of the 
orbit. In the other we will develop a similar transformation matrix via the 
three libration angles and the true position of the Moon as given by its right 
ascension and declination. These will then be equated and a solution obtained 
for the three libration angles. If we were interested only in determining the 
vector of a lunar landmark in the equatorial system as disturbed by the geo- 
metrical librations, the first transformation would be sufficient. However, 
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LUNAR ORBIT 


ECLIPTIC 


r MOON < i^EQUATOR 


Figure 3,- Orbital geometry of the mean Moon relative to the 
Earth and the equatorial inertial coordinate system. 


the physical librations are 
given, very conveniently, as a 
function of time, in terms of 
the three libration angles . 
Therefore, it is appropriate 
to compute the angles for the 
geometrical librations so that 
those determined for the 
physical librations may be 
added directly to them. The 
transformation matrix of the 
landmark is then formed from 
the total libration angles . 

The general geometry of 
the situation is shown in 
figure 3. The geocentric 
equatorial inertial coordinate 
system is shown as (IJK) , the 
obliquity of the ecliptic is 
e, the longitude of the lunar 
orbital node is fl, inclina- 
tion of the lunar orbit to the 


ecliptic is i, the argument 
of latitude is £ - Q. Note 
that i , the mean longitude of 
the Moon, is measured from the 
vernal equinox T along the 
ecliptic to the node and thence 



along the orbit to the position 
NORTH POLE of the mean Moon (not the true 

JNAR EQUATOR Moon) . The facts expressed by 

j Cassini's third law are shown 

in figure 4 for the instant 
the Moon occupies the ascending 
node. The inclination of the 
. lunar orbit relative to the 

ecliptic is i and that of 
% the lunar equator is I . The 

j poles of the three planes are 

sylk coplanar, with that of the 

ecliptic occupying the central 

\IN_UNAR ORBIT P° sition - The coplanarity of 
the poles is, of course, not 
ECLIPTIC disturbed by the progress of 
JT''''LUNAR EQUATOR the Moon in its orbit. 


NORTH POLE 
ECLIPTIC 


NORTH POLE 
LUNAR ORBIT 


Some time after the Moon 
has passed the node it is 
situated as shown in figure 3 


A close-up of the geometry at 


Figure 4.- Illustration of Cassini's Laws with 
the mean Moon at the orbital node. 
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NORTH POLE 

LUNAR AXIS this lunar position is shown 

in figure 5. The seleno- 
graphic coordinate system, 
which is fixed in the body of 
the Moon, is designated (ijk) . 
The i axis is located in the 
equatorial plane and is the 
origin of the longitude and 
latitude measurements, similar 
to the geographic system. The 
k axis is situated along the 
lunar rotational pole, and the 
j axis also lies in the plane 
of the lunar equator. The i 
axis is arbitrarily chosen as 
that radius of the mean Moon 
pointing to the center of the 
Earth when the mean Moon is at 
the node of the orbit. (Now, 
since the mean Moon and the 
true Moon coincide only at 
perigee and apogee, the i 
axis of the true Moon will 
point to the center of the 
Earth only when the true Moon 
is simultaneously located at the node, either ascending or descending, and at 
an apse, either perigee or apogee.) The radius vector between the Earth and 
the Moon (fig. 5) has rotated, in the orbital plane, from the node by the 
angle £ - Q. The i axis has also rotated from the nodal line by the angle 
£ - L but lies in the plane of the lunar equator. Our process now is to 
develop the transformations required to rotate the i axis into the radius 
vector to the Earth, as shown in figure 5, rotate this position by 180°, to 
point away from the Earth, and then, as shown in figure 3, to rotate along the 
orbital plane to the node, and then along the ecliptic to the vernal equinox. 
The transformations are given by the following equation: 



Figure 5.- Geometry between the mean lunar orbital 
equatorial planes and the terrestrial radius 
vector. 


x= [ex][^ z ][ix][(^)z][180 z ][(£-fi) z ] T [(I+i) x ][(£-fi) Jx' (Bl) 


The individual rotation matrices are here formed from the name of the angle 
and subscripted by the axis of rotation. The three rotation matrices, about 
the three axes, for the angle 0, are defined as: 
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The symbol x' represents any vector in the s elenographic coordinate system, 
while x is that same vector in the equatorial system. The first four 
matrices of equation (Bl) , in order from right to left, carry the i axis 
into the line of nodes (fig. 5), the equatorial plane into the orbital plane, 
the i axis into the radius vector, and then rotate it 180° so that it is now 
pointing directly away from the center of the Earth. The i axis is now 
alined with the radius vector from the Earth to the Moon (fig. 3), and the 
next four rotation matrices carry the i axis along the orbital plane to the 
node, rotate the orbital plane into the plane of the ecliptic, the i axis 
along the ecliptic to the vernal equinox, and the ecliptic into the plane of 
the Earth's equator. 

The equation may now be simplified by factoring the rotation matrix 
[ 1 80 z ] out to the left. 1 In so doing, all rotations about the z axis are 
unaffected and any rotations about the x or y axis are transposed. The 
result and the further steps in the simplifications are given below. 

x = [180 z ][ ex ] T [fi z ][i x ] T [(iL-n) z ][(£-S2) z ] T [(I+i) x ][(^-fi) z ]x''| 

x = [180 z ][c x ] T [fi z ][i x ] T [i x ][I x ][(£-Q) z ]x' > (B2) 

x= [180 z ][e x ] T [Q z ] [I x ][U-«) z ]x' f 


1 This is a legal operation, because of the special nature of the [180 z ] 
matrix. 
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Equation (B2) transforms any vector in the s elenographic coordinate system 
into a vector in the equatorial system. As mentioned above, if we were 
interested only in the effects of the geometrical librations this equation 
would suffice. 


MOON 

NORTH POLE 



TO EARTH 
CENTER 


/EQUATOR 


CELESTIAL 
NORTH POLE 



-T 


EARTH EQUATOR 


'"'"•PLANE OF THE 
CELESTIAL EQUATOR 


(b) The position of the Moon in the geocentric equatorial 

system. 

Figure 6.- The true Moon as it appears in the heavens to an 
observer on the Earth. 


We will now develop the 
transformation matrix via the 
libration angles and the true 
position of the Moon as given 
by the right ascension and 
declination. In figure 6, the 
Moon is shown as it appears in 
the heavens to an observer on 
Earth. An hour circle of the 
equatorial coordinate system 
is shown passing through the 
apparent center C of the 
Moon, and is marked N at its 
northern part and S at its 
southern part. A small circle, 
known as a declination circle 
is marked E for its eastern 
part and W for its western 
part. The central meridian, 
in the s elenographic coordin- 
ate system, is shown passing 
through C and the north pole 
of the Moon. The position 
angle of this meridian is 
designated P 0 , and normal to 
it is the lunar equator. The 
point 0 is the origin of 
s elenographic coordinates and 
extending from it is the 
principal axis i. To an 
observer on the Moon at the 
point C, the center of the 
Earth appears to be in the 
zenith. The selenographic 
coordinates of the point C 
are given as longitude X and 
latitude g, both positive as 
shown. The three angles A, 
g, and P Q are the three libra- 
tion angles in longitude, lati- 
tude, and position angle. In 
addition, the position of the 
center of the Moon is given as 
a in right ascension and 6 
in declination relative to the 
equator and equinox of our 
fundamental system. We are 
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now in a position to formulate our transformation matrix, which we may write 
as : 


x= [a z ][S y ][180 z ][P Ox ][3 y ] T [X z ] T x' (B3) 

Here, reading the rotation matrices from right to left, we have rotated the 
i axis through the longitude angle A and the latitude angle 8 until it is 
pointing in the direction of the center of the Earth, and then rotated the 
plane of the lunar equator through the angle P Q until it is coincident with 
the declination circle EW and parallel to the equator. Next we have rotated 
the i axis until it is pointing directly away from the Earth, and then 
brought it into coincidence with the vernal equinox through the successive 
angles 6 and a. As with the previous development we may factor the matrix 
[180 z ] out to the left, producing 

x = [180 z ][a z ][6 y ] T [P Ox ][8 y ] T [A z ] T x' (B4) 

We now have two different expressions (eqs. (B2) and (B4)) containing 
the same transformation matrix which we may equate and solve for the unknown 
libration angles. 

[Po x ][fV] T U z ] T = [6 y ][a z ] T [£ x ][U z ][I x ][0l - n) z ] (B5) 

The quantities on the left are, of course, our unknowns and all the quantities 
on the right are available for computation. The only problem left is that of 
unscrambling the three angles from the left side of the expression. This may 
be readily accomplished by setting 

[A] = [P 0x ][8 y ] T [A z ] T (B6) 

and then expanding the right-hand side to provide: 


[A] = 


(cos g cos X) 

(cos P 0 sin X + sin P 0 sin B cos X) 
(-sin P 0 sin X - cos P 0 sin g cos X) 


(cos 6 sin X) 

(cos P 0 cos X + sin P 0 sin g sin X) 
(sin P 0 cos X - cos P Q sin B sin X) 


(sin B) 

(-sin P Q cos B) 
(cos P 0 cos B) 


(B7) 


The matrix [A] may be computed from the right-hand side of equation (B5) and 
we may then find our three libration angles from the elements of this matrix. 


tgA 
sin 6 

tgPo 


A 1 2 
An" 

A 1 3 
A 


23 


v 33 


(B8) 
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The values of the right ascension a and the declination 6 of the 
Moon have been precomputed and are available in tabular form or on magnetic 
tape suitable for an electronic computer. The other quantities are either 
constant or have been expressed as functions of time. These quantities 
I, e, l, and S3 have been taken from reference 9 and are given below: 

I = 1°32 ! 1 

£ = 23°27 ' 08726 - 467845T - 0V0059T 2 + 0V00181T 3 Epoch = 1900.0 

= 1900 Jan 0.313 

l = 270?434164 + 13?1763965268d - 0?0000850D 2 Epoch = 1900 Jan 0^5 ET 

= JD 2415020.0 

S3 = 259?183275 - 0?0529539222d + 0?0001557D 2 
where 

d JD - 2415020.0 
D d/ 10 , 000 

T (d + 0.187)/36525.0 

JD Julian Day number 

The lunar geometrical librations may be computed from the above simple 
equations, which obviates the need to store the values on tape and then 
interpolate them to the desired time. 

The physical librations may be computed from relatively simple equations 
that provide the increments AA, Ag, and AP 0 to be added to the libration 
angles discussed above. From reference 9, 

AA = 0?003 sin(£~r ' ) - 0?005 sin 2 (I” -S3) - 0?016 sin g Q + AP 0 sin g\ 

Ag = M cos A + N sin A > (B9) 

AP Q = (M sin A - N cos A)sec g / 

where 

M = 0?040 sin(T ' - S3) - 0?003 sin(£ - S3) 

N = 0?020 cos(r' - S3) + 0?003 cos (£ - S3) 

T' = 334?329556 + 0?1114040803d - 0?0007739D 2 Epoch = 1900 Jan 0?5 ET 

= JD 2415020.0 

g = 358?475845 + 0?985600267d - 0?0000112D 2 

s o 
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Here £ and Q are as defined above, I” is the mean longitude of the 
lunar perigee, and g Q is the mean anomaly of the Sun; M and N are merely 
auxiliary quantities. 

The total libration angles are then computed by adding the increments 
found from equations (B9) above to the geometrical libration angles computed 
from equations (B8) . 


X ’ = X + AX 
S' = 3 +■ A3 

V = P o + AP o 

and the vector of the landmark in equatorial coordinates is then found by 
substituting into the equation (B4) , which becomes 

x = [180 z ][a 2 ][6 y ] T [P^][B^] T [X'] T x' (BIO) 

It will be noted that the equations for e, SL , fi, F', and g are 
referenced to a particular epoch and that they are given in terms of the time 
variables d, D, and T. The quantity d is the number of Julian days that 
have elapsed since the epoch 1900 January 0^5 (or Julian day number 2415020.0). 
In Civil time, this epoch is December 31, 1899, at Greenwich noon. The 
quantity D is the same as d but expressed in terms of 10,000 Julian 
days. Similarly, T is essentially the same as d expressed in terms of 
Julian centuries of 36525 days. The quantity 0.187, in the defining equation 
for T, is an adjustment for the difference in epochs between the one equation 
containing T (the equation for e) , whose epoch is 1900 January 0^313, and 
the other equations containing d and D whose epoch is 1900 January 0^5. 

The difference between ephemeris and universal time at these fundamental 
epochs, which is -3^9, has been considered to be sufficiently small to be 
ignored for the present purposes. It should be carefully noted that all of 
these equations provide values of their respective variables relative to the 
mean equator and equinox of date. As a consequence our landmark vector found 
from equation (B10) is represented in this same coordinate system and must be 
precessed to the fundamental system before further use is made of it. 
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APPENDIX C 


DEVELOPMENT OF EQUATIONS FOR WINDOW REFRACTION 


An observer in a presurrized spacecraft may be required take his 
navigation observations through a window, or a set of window panes. The 
windows will undoubtedly be deflected by the differential pressure existing 
between the interior and the exterior, and any light ray traversing the win- 
dow will be refracted due to both the difference in pressure and the deflec- 
tion of the window. The external ray to the object sighted may always be 
considered to be approximately known, since we assume that a navigation sight 
to a known object is being made. 

The first step to an understanding of the process is to realize that an 
entire bundle of parallel rays from the object intercepts the window at all 
points on the surface. And that, even though all these rays are transmitted 
through the window, only one ray enters the eye of the observer, or the 
instrument he may be using. The real problem, then, is that of discovering 
this particular ray. Since the mathematical description of optical phenomenon 
is reversible we will consider a ray to be emerging from the eye of the 
observer, extending it to the window, refracting it at the window, and then 
continuing it on to the body under observation. When a higher order function 
is used to represent the surface of the refracting window, the problem is 
mathematically unsolvable in an analytical form and recourse must be had to 
the process of approximating the solution by some method of iteration. The 
method that has been set up to be used at Ames in a digital computer program 
(see ref. 18) is that of first assuming the direction of a ray from the eye 
of the observer toward the window. This ray is then carried forward by an 
iterative process until it intersects the surface of discontinuity, at which 
point the refracted ray beyond the surface is determined. The process is 
repeated in sequence for each surface of the window. When the last surface 
has been considered, the direction of the refracted ray beyond the window is 
then known and can be compared with the direction of the body specified as 
having been observed. The rotation vector between the refracted ray and the 
true ray to the object is then determined, modified by the amount of the 
refraction and applied to the original, assumed, interior ray. The whole 
process is then repeated as many times as necessary until the refracted ray 
and the true ray coincide. At this point the last modified ray is the one 
desired. As many steps in this loop may be added as desired to account for 
multiple refracting surfaces. 

It is necessary then, to formulate equations that will represent: the 

three-dimensional window surface, the normal to the surface, the refraction 
of the ray, and the iterations that may be required. All are to be repre- 
sented in vector-matrix form for ease in computing these three-dimensional 
entities . 

The computation of the refraction of the rays of an object observed 
through the windows of a pressurized space capsule is complicated by the fact 
that the window shape must necessarily be dictated by considerations of space 
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available for windows and their placement to accommodate critical modes of 
flight and maneuver. As a consequence, the windows, as bowed under a differ- 
ential pressure, are almost certain to have a distinct lack of symmetry, and 
the rays to be considered are almost certain to be nonradial; previously 
developed ray trace analyses therefore are not applicable. 


The purpose of this appendix is to develop an applicable ray-trace 


analysis. The development 
coordinates, since that is 


A 

K 


T" 



Figure 7.- Orientation of the spacecraft coordinate 
system (I'J'K 1 ) with respect to the equatorial 
inertial coordinate system (IJK) . 


A ' 

K 



Figure 8.- Orientation of the window coordinate system (ijk) 
with respect to the spacecraft coordinate system (I'J'K'). 


of this analysis will be carried out in window 
the coordinate system in which the window is most 

easily described and is thus 
most convenient to use when 
iterations are required. 

One vector, however, must be 
transformed into the window 
coordinate system before the 
process may be carried 
through. This is the unit 
vector of the observed celes- 
tial body, which is provided 
in the equatorial inertial 
coordinate system. 

The geometry of the 
transformation may be under- 
stood by referring to fig- 
ures 7 and 8. In figure 7, 
the fundamental inertial 
coordinate^s^stem is desig- 
nated as (IJK) , with the 
principal axis indicated by 
the direction to the vernal 
equinox T. The body axis 
of the space capsule is 
given as (I'J'K'), and the 
orientation angles relative 
to the fundamental system as 
<j>, ip, and 0. The angles 4> 
and ij) are similar to right 
ascension and declination, 
and they describe the inter- 
section of the principal 
body axis with the celestial 
sphere, while 0 is the 
roll angle of the body. Fig- 
ure 8 shows the orientation 
of the window coordinate 
system (ijk) with respect to 
the body axis system (I'J'K*). 
The space capsule is shown 
as conical in shape, since 
that seems to be the type of 
configuration being consid- 
ered for space capsules of 
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the near future. The half-cone angle of the body is shown as C, and the 
roll angle of the window out of the reference plane of the space capsule is 
given as W. From figure 8, it can be seen that the transformation from the 
body axis to the window axis system may be written as 



where 


[T WB ] = [c y ][w x ] T 

And the transformation from the equatorial system to the body axis system may 
be written as 


(i’l 


fi) 

1 



r 

’ = [Tbe] ' 

d 

UJ 


k/ 


where 

[Tee! - [e x ] T [* y ] T [f z ] T 

Then the complete transformation from the equatorial system to the window 
coordinate system may be given as 


n 

n 

j / = [ t Wb] [ T Be] ' 

j 

U 

,K> 


Or in terms of a unit vector Sg in the equatorial system and a unit vector 
% in the window coordinate system 

S W = [C y ] [W x ] T [0 x ] T [i|] y ] T [(j) z ] T S E (Cl) 

The notation used here for unit angular rotation matrices is the same as that 
given in appendix B. 


Next, one must develop a vector-matrix equation representing the three- 
dimensional refracting surface which may be accomplished in two steps. The 
first of these requires that the general shape of the surface be generated by 
means of a structural deformation-type program that makes use of the physical 
characteristics of the window, the shape and thickness of the window, and the 
pressure differential across the window. Programs of this type are available 
for use; one of these is described in reference 19. The output of this 
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program is the deflection of the window at discrete points over its surface 
for any pressure differential. Although this type of information represents 
the surface with sufficient accuracy, it is in a nonanalytical form, and if it 
is used directly a great deal of interpolation is required since the rays will 
invariably fall between the points for which the discrete data are available 
Consequently, as a second step, it is considered desirable to convert the 
data to an analytical form that is smooth and tractable over the portion of 
the window that will come under consideration, and in addition, may be 
expressed in the desired vector-matrix form. This may be accomplished by 
fitting polynomials to cuts through the surface, and fitting the coefficients 
of these equations with polynomials in the remaining variable. The edges of 
the windows, as installed in the space capsule, may have clamped or unclamped 
edges. Our equation should be capable of representing the highest order 
curve thus generated. This would be the case with clamped edges, which pro- 
duce curves with two points of inflection. A fourth-degree equation, used 
here, will represent this type of surface. This technique can be used for 
higher order equations, for which the necessary modifications will be 
apparent . 


Let us now consider that we have taken cuts through our surface in 
planes parallel to the yz plane at constant values of the x coordinate. 
Fitting the resultant deflections, in the z coordinate, by least squares we 
now have several equations of the form 

z = A x y 4 + B x y 3 + C x y 2 + D x y + E x (C2) 

where the coefficients are functions of x. Equation (C2) may now be factored 
into the form 


z = (y 4 y 3 y 2 y 1) 



(C3) 


Since the coefficients are functions of x and are known from the above 
mentioned least-squares fits to the original data, we may represent them as 


A x 

= m n x 4 

+ 

mi 2 x 3 

+ 

mi 3 x 2 

+ 

mj 4 x 

+ 

mi 5 ' 

B X 

= m 2 iX 4 

+ 

m 2 2 x 3 

+ 

m 23 x 2 

+ 

m 24 x 


m 2 5 

c x 

= m 31 x 4 

+ 

m 32 x 3 

+ 

m 33 x 2 

+ 

m 34 x 

+ 

m 3 5 > 

D X 

= m 41 x 4 

+ 

m 42 x 3 

+ 

m 43 x 2 

+ 

m 44 x 

+ 

m 4 5 

E x 

= m 51 x4 

+ 

m 52 x 3 

+ 

m 53 x2 

+ 

m 54 x 

+ 

m 5 5 ^ 


(C4) 


The values on the left of these equations are known, as is the variable x; 
hence we may fit these coefficients with a least-square process also. 
Equations (C4) may now be factored into the form 


37 




Combining equations (C3) and (C5) yields 

z = Y 7 [M]X 

where 


[M] 


mn 

m l 2 

m l 3 

m 14 

m 21 

m 2 2 

m 2 3 

m 2 4 

m 31 

m 32 

m 33 

m 34 

m 41 

m 4 2 

m 4 3 

m 44 

m 51 

L. 

m 52 

m 53 

m 54 



(C6) 


(C7) 


X T = (x 4 x 3 x 2 x 1) 

Y T = (y 4 y 3 y 2 y 1) 

Equation (C7) is known as the surface matrix since it represents the shape of 
the surface of the bowed window under a differential pressure of unity. When 
the differential pressure is other than unity, equation (C6) is multiplied 
by the gage pressure Ap, since we may consider the deflection to be a linear 
function of the pressure. 

z = ApY T [M]X (C8) 

Although the quantities Y and X are not true geometrical vectors, since 
their elements are not independent, the form of the equation suits our pur- 
poses because an expansion of equation (C6) or (C8) produces the mixed 
polynomial that we require. 
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As indicated in the above development, the [M] matrix would be developed 
by first fitting curves as a function of y by least squares at several con- 
stant values of x, and then the coefficients from these first fits would 
themselves be fitted by least squares as a function of x. A more expeditious 
method may be that of forming matrices from the various given values of x, y, 
and z, and solving for the [M] matrix in a single matrix least-squares solu- 
tion. This method may be readily developed as follows from equation (C6) . 

z = Y T [M]X 


By stacking several of these equations under each other for a constant value 
of x and factoring into vector-matrix form we have 

Z = [Y] [M]X 

Now if we stack several of these equations alongside of each other at constant 
values of y, again factoring into the vector-matrix form, we obtain 

[Z] = [Y] [M] [X] 


The standard method of forming the normal equations in matrix notation may now 
be applied by multiplying through the equation from the left by [Y] T and from 
the right by [X]' r . Thus 

[Y] T [Y][M][X][X] T = [Y] T [Z][X] T 


Now if sufficient, redundant, independent observations are available for the 
least-squares reduction, the quantities [Y]^[Y] and [X] [X] 1 are square, 
nonsingular, and symmetrical, but nonorthogonal . Thus if we solve for the 
matrix [M] we have 


[M] = ([Y] T [Y])- 1 [Y] T [Z][X] T ([X][X] T )‘ 1 


CC9) 



Figure 9.- Geometry of the scheme for 
iterating from the base plane to the 
surface of discontinuity. 


This is the equation used to deter- 
mine the surface matrix [M] when the 
description of the surface is 
supplied in a rectangular tabular 
form. 

Now, having the mathematical 
description of the window surface, we 
can determine the intersection of the 
assumed interior ray with this sur- 
face. This is an iterative procedure, 
a and is illustrated in figure 9. Here 
the window coordinate system is given 
as (ijk), the vector from the origin 
locating either the eye of the 
observer or some reference of the 
instrument used, is given as S, some 
initial ray vector as Vi, which ter- 
minates in the xy plane, and a vec- 
tor in the xy plane closing with 
Vj as d. Let 
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S = S x i + Syj + S z k 
Vi = ai + 3j + yk 


Vi = ViV 


where V! is unknown. Now d is of the form 


and from figure 9 


decomposing 


then from the third equation 


and from the first two 


d = xji + + Ok 


d = S + ViV 


X 1 = S x + V x a 


Yl = Sy + 

0 = S z + V X Y 


Vi = - 


X 1 ~ S X " y 

n - s y - 7 s z 


ceil) 


Now, if we set 


AVi = 0 


CC 12) 


compute 


Zl = ApY: [M] X ! 


(C13) 


and project the increment of Zl above the Vj vector onto the V direction 


Rl = ( Z i - AV X V • k)k • V = ( Zl - AViy)y 


(C14) 


Now form 


AV2 — AVj + Rj 


CC15) 


and compute 
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(C16) 


^ ^ V 

X 2 = xi + R^V • i = Xx + R^a 
y 2 = yi + RjV • j = yi + 

We may now loop back to the equation for z (eq. (C13)) and continue to 
compute until |r| < e. At this time, the last values computed for x, y, 
and z are the coordinates of the point of intersection of our ray with the 
refracting surface, to the accuracy which we have specified for e. The basic 
scheme of this iteration is that by successively projecting the value of 
(z - RV • k) onto the vector V we approach the point of intersection with 
the window surface. This system appears to be stable for all continuous 
smooth surface functions . 

When the coordinates of the point of intersection of the ray vector with 
the surface has been established, the unit vector normal to the surface (N) 
may be determined from the gradient of the scalar function F 


N 


VF 

VFj 


(C17) 


where 


F(x,y,z) = z - ApY [M]X = 0 



+ 



k 



(4x 3 3x 2 2x 1 0) 



(4y 3 3y 2 2y 1 0) 
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and 


VF 



After the unit vector normal to the surface has been established, the 
refracted ray may be computed by an equation which will be developed with the 
help of figure 10. The unit vector N is the normal to the surface of 

discontinuity, on either side of which exist 
a medium of transmission whose indices of 
refraction, relative to free space, are 
n 2 and n^ . The incident ray has unit vector 
I and angle of incidence i, while the 
refracted ray has unit vector R and refrac- 
tion angle r. The two principal precepts of 
geometrical optics are that the vectors 
I, R, and N are coplanar, and that Snell's 
Law holds : 



Ik 



sm 1 


"2 

n l 


sin r 

From the requirement of coplanarity, we may 
write 


I x N = sin i c 


and 


Figure 10.- Refraction at a surface of 
discontinuity, 


R x N = sin r c 

where c is a unit vector normal to i, R, 
and N. The elimination of c from the two 
equations leaves 


(R x N) 


sin r 
sin i 


(I x N) 


n 2 


(I x N) 


If now, we cross through this equation with N from the left, we may solve 
for the vector R, which gives us 


R = 





[1 - (I • N) 2 ] 



(I ■ N) 


N 


(CIS) 


Although this equation has been developed on the basis of refraction, and we 
will use it for that purpose, it applies equally well to reflection by merely 
setting ni and n 2 equal to each other and taking the negative sign on the 
radical . 
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The star unit vector, already transformed into the window coordinate 
system, is now compared with the refracted ray R. (For simplicity, we will 
call the star vector in the window coordinate system S.) This may be 
accomplished by comparing the angle with an arbitrary assigned limit, such as 

|R x s| <_ e (C19) 

Now if the refracted ray R is not alined with the star vector (fig. 11), it 
is necessary to modify the internal vector V and continue with another 

iteration. The modification 
of the old internal vector 
may be accomplished in the 
following manner. Two 
assumptions are made. The 
first is that the internal 
vectors and the external vec- 
tors are coplanar; that is 

R x S _ ^old x ^new 
|R x s| |V old x v new 

The second is that the magni- 
tude of the angle change is 
modified by the ratio of the 
change in the angle of inci- 
dence to the change in the 
Figure 11.- Rotation of R into the star vector and generation angle of refraction. If G 
°f Vnew f rom v oid- represents this ratio, then 

G|R x s| = [V old x V new | 

When the second of these equations is substituted into the first the result 
is 

x V n ew = G(R x §) 

This may be solved for V new by crossing through the equation from the 
right with V 0 ]_ d . The result is 

Vew = °[(R x S) x V old J + V Qld 

with the assumption that 




v new ' v old " 1 

Since this assumption will destroy the unit magnitude of the vector, the 
unit vector may be found from 


G[(R x S) x V Qld ] f V Qld 
G[(R x S) x V old ] + V Qld 


(C20) 
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The gain factor G may be derived simply by starting with the equation 


sin i _ n 2 
sin r nj 


and differentiating, we have 

ri 2 

cos i di = — cos r dr 
n l 


and this we may reduce to 


/di\ = /W\ (R • N) 

y dr y ( i . n) 


which is the gain factor through a single refracting surface. Now if several 
refracting surfaces are in the path of the ray, then we may write, with 
sufficient approximation 



(C21) 


since the change in the incident ray at a surface is equal to the change in 
the refracted ray at the preceding surface . 

The entire process in its proper order may be conveniently listed in 
step form as follows: 

1. Transform the observed star unit vector from equatorial inertial 
coordinates to window coordinates using equation (Cl) ■ 

2. Assume any convenient value of the internal unit vector V (values 

of a, 3, and y in eq . (CIO)). A convenient value is a = 3 = 0, and y = 1. 

3. Compute the intersection of this vector with the plane of the 
undeflected window using equation (Cll). 

4. Set Vi = 0, equation (C12). 

5. Using equations (C13) to (C16) , iterate to find the intersection of 
the ray with the deflected surface of the window. 

6. Compute the normal to the surface using equation (C17). 

7. Compute the refracted ray using equation (C18) . 

8. Repeat steps 3 to 7 for each of the refracting surfaces. 

9. Using equation (C19) , determine the angle between the free-space 
refracted ray and the star vector determined in step 1. 
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10. If the angle determined in step 9 is too large, compute a new 
internal unit vector using equation (C20) , including the gain factor computed 
from equation (C21) . 

11. Iterate the entire process from steps 3 through 10 until the angle 
determined in step 8 is small enough. 

12. The final refracted ray is then coincident with the star ray and 
the proper initial internal ray vector has been determined. 

This theory has been programmed in FORTRAN IV for operation on an 
IBM 7094 digital computer. The operation of the program is completely dis- 
cussed in reference 14. An indication of the accuracy of the theory, as 
programmed, relative to thin lens theory may be gained by comparing the fig- 
ures in table II. These were obtained by computing the deflection of 


TABLE II.- COMPARISON BETWEEN DEFLECTION ANGLE COMPUTED BY THIN 
LENS THEORY AND THE PRESENT THEORY AND COMPUTER PROGRAM 


Angle of 
incident 
of ray 

Vertical 
location 
of ray 

Deflection angle 
from thin 
lens theory 
(sec of arc) 

Deflection angle 
from present 
theory and program 
(sec of arc) 


2.5 

12.176 

12.182 


1.25 

6.088 

6.094 

0 

0 

0 

.006 


-1.25 

-6.088 

-6.082 


-2.5 

-12.176 

-12.170 


par.a-axial rays passing through a thin spherical lens 6 inches in diameter, 
0.38 inch thick at the edges with two spherical surfaces of 11.60xl0 - ' S inch 
sagittal height. The dimensions of this lens corresponds to a portion of a 
simulated spacecraft window for which other test data were computed. The 
equation, in thin lens theory, used as a check was 



where 
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n 2 / n l = 1-458 

D = 0.380 in. 

R = 38,793.10 in. (radius of curvature) 

The matrix [M] used in equation (C7) for the computation in the present theory 
was 


64 


JL_ D“ 5 
16 R 


1 R-3 
8 R 


0 


[M] = 


3 r>“5 


16 

0 


R D 0 


jL D ” 3 
4 R 


i R_1 


I n-3 
8 R 


I 0-1 


This was obtained from the equation for a spherical surface of sagittal 
height h 

z = h - R + [R 2 - x 2 - y 2 ] 

by expanding in the binomial theorem. It is apparent from the two right-hand 
columns of table II that the two methods of computation differ by a maximum 
of only 0'.'006. This difference could well be due to round off errors in the 
computation. 
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